{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Using `astropy.coordinates` to Match Catalogs and Plan Observations\n",
    "\n",
    "## Authors\n",
    "Erik Tollerud, Kelle Cruz\n",
    "\n",
    "## Learning Goals\n",
    "* TBD\n",
    "\n",
    "## Keywords\n",
    "coordinates, catalog matching, observational astronomy, astroquery\n",
    "\n",
    "## Summary\n",
    "In this tutorial, we will explore how the `astropy.coordinates` package and related astropy functionality can be used to help in planning observations or other exercises focused on large coordinate catalogs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Python standard-library\n",
    "from urllib.parse import urlencode\n",
    "from urllib.request import urlretrieve\n",
    "\n",
    "# Third-party dependencies\n",
    "from astropy import units as u\n",
    "from astropy.coordinates import SkyCoord\n",
    "from astropy.table import Table\n",
    "import numpy as np\n",
    "from IPython.display import Image\n",
    "\n",
    "# Set up matplotlib and use a nicer set of plot parameters\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Section 0: Describing on-sky locations with `coordinates`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's start by considering a field around the picturesque Hickson Compact Group 7.  To do anything with this, we need to get an object that represents the coordinates of the center of this group."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center = SkyCoord(9.81625*u.deg, 0.88806*u.deg, frame='icrs')\n",
    "hcg7_center"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-info\">\n",
    "\n",
    "**Note:** If you already worked through [Coords 1: Getting Started with Coordinates](http://learn.astropy.org/rst-tutorials/Coordinates-Intro.html), feel free to skip to [Section 1](#Section-1:).\n",
    "\n",
    "</div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In Astropy, the most common object you'll work with for coordinates is `SkyCoord`.  A `SkyCoord` can be created most easily directly from angles as shown below.  It's also wise to explicitly specify the frame your coordinates are in, although this is not strictly necessary because the default is ICRS. \n",
    "\n",
    "(If you're not sure what ICRS is, it's basically safe to think of it as an approximation to an equatorial system at the J2000 equinox)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "SkyCoord will also accept string-formatted coordinates either as separate strings for ra/dec or a single string.  You'll have to give units, though, if they aren't part of the string itself."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "SkyCoord('0h39m15.9s', '0d53m17.016s', frame='icrs')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "SkyCoord('0:39:15.9 0:53:17.016', unit=(u.hour, u.deg), frame='icrs')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If the object you're interested in is in [SESAME](http://cdsweb.u-strasbg.fr/cgi-bin/Sesame), you can also look it up directly from its name using the `SkyCoord.from_name()` class method<sup>1</sup>. Note that this requires an internet connection.  It's safe to skip if you don't have one, because we defined it above explicitly.\n",
    "\n",
    "If you don't know what a class method is, think of it like an alternative constructor for a `SkyCoord` object -- calling `SkyCoord.from_name()` with a name gives you a new `SkyCoord` object. For more detailed background on what class methods are and when they're useful, see [this page](https://julien.danjou.info/blog/2013/guide-python-static-class-abstract-methods)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center = SkyCoord.from_name('HCG 7')\n",
    "hcg7_center"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This object we just created has various useful ways of accessing the information contained within it.  In particular, the ``ra`` and ``dec`` attributes are specialized [Quantity](http://docs.astropy.org/en/stable/units/index.html) objects (actually, a subclass called [Angle](http://docs.astropy.org/en/stable/api/astropy.coordinates.Angle.html), which in turn is subclassed by [Latitude](http://docs.astropy.org/en/stable/api/astropy.coordinates.Latitude.html) and [Longitude](http://docs.astropy.org/en/stable/api/astropy.coordinates.Longitude.html)).  These objects store angles and provide pretty representations of those angles, as well as some useful attributes to quickly convert to common angle units:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "type(hcg7_center.ra), type(hcg7_center.dec)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center.dec"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center.ra"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center.ra.hour"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now that we have a `SkyCoord` object, we can try to use it to access data from the [Sloan Digitial Sky Survey](http://www.sdss.org/) (SDSS).  Let's start by trying to get a picture using the SDSS image cutout service to make sure HCG7 is in the SDSS footprint and has good image quality.\n",
    "\n",
    "This requires an internet connection, but if it fails, don't worry: the file is included in the repository so you can just let it use the local file``'HCG7_SDSS_cutout.jpg'``, defined at the top of the cell.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "impix = 1024\n",
    "imsize = 12*u.arcmin\n",
    "cutoutbaseurl = 'http://skyservice.pha.jhu.edu/DR12/ImgCutout/getjpeg.aspx'\n",
    "query_string = urlencode(dict(ra=hcg7_center.ra.deg, \n",
    "                              dec=hcg7_center.dec.deg, \n",
    "                              width=impix, height=impix, \n",
    "                              scale=imsize.to(u.arcsec).value/impix))\n",
    "url = cutoutbaseurl + '?' + query_string\n",
    "\n",
    "# this downloads the image to your disk\n",
    "urlretrieve(url, 'HCG7_SDSS_cutout.jpg')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now lets take a look at the image."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Image('HCG7_SDSS_cutout.jpg')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Section 1: Using `coordinates` and `table` to match and compare catalogs"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "At the end of the last section, we determined that HCG7 is in the SDSS imaging survey, so that means we can use the cells below to download catalogs of objects directly from the SDSS. Later on, we will match this catalog to another catalog covering the same field, allowing us to make plots using the combination of the two catalogs."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will access the SDSS SQL database using the [astroquery](https://astroquery.readthedocs.org) affiliated package.  This will require an internet connection and a working install of astroquery. If you don't have these you can just skip down two cells, because the data files are provided with the repository. Depending on your version of astroquery it might also issue a warning, which you should be able to safely ignore."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from astroquery.sdss import SDSS\n",
    "sdss = SDSS.query_region(coordinates=hcg7_center, radius=20*u.arcmin, \n",
    "                         spectro=True, \n",
    "                         photoobj_fields=['ra','dec','u','g','r','i','z'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`astroquery` queries gives us back an [astropy.table.Table](http://docs.astropy.org/en/stable/table/index.html) object. We could just work with this directly without saving anything to disk if we wanted to. But here we will use the capability to  write to disk. That way, if you quit the session and come back later, you don't have to run the query a second time.\n",
    "\n",
    "(Note that this won't work fail if you skipped the last step.  Don't worry, you can just skip to the next cell with ``Table.read`` and use the copy of this table included in the tutorial.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sdss.write('HCG7_SDSS_photo.dat', format='ascii')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If you don't have internet, you can read the table into python by running the cell below.  But if you did the astroquery step above, you could skip this, as the table is already in memory as the `sdss` variable."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sdss = Table.read('HCG7_SDSS_photo.dat', format='ascii')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ok, so we have a catalog of objects we got from the SDSS.  Now lets say you have your own catalog of objects in the same field that you want to match to this SDSS catalog.  In this case, we will use a catalog extracted from the [2MASS](http://www.ipac.caltech.edu/2mass/).  We first load up this catalog into python."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "twomass = Table.read('HCG7_2MASS.tbl', format='ascii')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now to do matching we need `SkyCoord` objects.  We'll have to build these from the tables we loaded, but it turns out that's pretty straightforward: we grab the RA and dec columns from the table and provide them to the `SkyCoord` constructor.  Lets first have a look at the tables to see just what everything is that's in them."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sdss # just to see an example of the format"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "twomass # just to see an example of the format"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "OK, looks like they both have ``ra`` and ``dec`` columns, so we should be able to use that to make `SkyCoord`s.\n",
    "\n",
    "You might first think you need to create a separate `SkyCoord` for *every* row in the table, given that up until now all `SkyCoord`s we made were for just a single point.  You could do this, but it will make your code much slower.  Instead, `SkyCoord` supports *arrays* of coordinate values - you just pass in array-like inputs (array `Quantity`s, lists of strings, `Table` columns, etc.), and `SkyCoord` will happily do all of its operations element-wise."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "coo_sdss = SkyCoord(sdss['ra']*u.deg, sdss['dec']*u.deg)\n",
    "coo_twomass = SkyCoord(twomass['ra'], twomass['dec'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note a subtle difference here: you had to give units for SDSS but *not* for 2MASS.  This is because the 2MASS table has units associated with the columns, while the SDSS table does not (so you have to put them in manually).\n",
    "\n",
    "Now we simply use the ``SkyCoord.match_to_catalog_sky`` method to match the two catalogs. Note that order matters: we're matching 2MASS to SDSS because there are many *more* entires in the SDSS, so it seems likely that most 2MASS objects are in SDSS (but not vice versa)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "idx_sdss, d2d_sdss, d3d_sdss = coo_twomass.match_to_catalog_sky(coo_sdss)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "``idx`` are the indecies into ``coo_sdss`` that get the closest matches, while ``d2d`` and ``d3d`` are the on-sky and real-space distances between the matches. In our case ``d3d`` can be ignored because we didn't give a line-of-sight distance, so its value is not particularly useful.   But ``d2d`` provides a good diagnosis of whether we actually have real matches:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.hist(d2d_sdss.arcsec, histtype='step', range=(0,2))\n",
    "plt.xlabel('separation [arcsec]')\n",
    "plt.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ok, they're all within an arcsecond that's promising.  But are we sure it's not just that *anything* has matches within an arcescond?  Lets check by comparing to a set of *random* points.\n",
    "\n",
    "We first create a set of uniformly random points (with size matching `coo_twomass`) that cover the same range of RA/Decs that are in `coo_sdss`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ras_sim = np.random.rand(len(coo_twomass))*coo_sdss.ra.ptp() + coo_sdss.ra.min()\n",
    "decs_sim = np.random.rand(len(coo_twomass))*coo_sdss.dec.ptp() + coo_sdss.dec.min()\n",
    "ras_sim, decs_sim"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we create a `SkyCoord` from these points and match it to `coo_sdss` just like we did above for 2MASS.\n",
    "\n",
    "Note that we do not need to explicitly specify units for `ras_sim` and `decs_sim`, because they already are unitful `Angle` objects because they were created from `coo_sdss.ra`/`coo_sdss.dec`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "coo_simulated = SkyCoord(ras_sim, decs_sim)  \n",
    "idx_sim, d2d_sim, d3d_sim = coo_simulated.match_to_catalog_sky(coo_sdss)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now lets plot up the histogram of separations from our simulated catalog so we can compare to the above results from the *real* catalog."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.hist(d2d_sim.arcsec, bins='auto', histtype='step', label='Simulated', linestyle='dashed')\n",
    "plt.hist(d2d_sdss.arcsec, bins='auto', histtype='step', label='2MASS')\n",
    "plt.xlabel('separation [arcsec]')\n",
    "plt.legend(loc=0)\n",
    "plt.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Alright, great - looks like randomly placed sources should be more like an arc*minute* away, so we can probably trust that our earlier matches which were within an arc*second* are valid.  So with that in mind, we can start computing things like colors that combine the SDSS and 2MASS photometry."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rmag = sdss['r'][idx_sdss]\n",
    "grcolor = sdss['g'][idx_sdss] - rmag\n",
    "rKcolor = rmag - twomass['k_m_ext']\n",
    "\n",
    "plt.subplot(1, 2, 1)\n",
    "plt.scatter(rKcolor, rmag)\n",
    "plt.xlabel('r-K')\n",
    "plt.ylabel('r')\n",
    "plt.xlim(2.5, 4)\n",
    "plt.ylim(18, 12) #mags go backwards!\n",
    "\n",
    "plt.subplot(1, 2, 2)\n",
    "plt.scatter(rKcolor, rmag)\n",
    "plt.xlabel('r-K')\n",
    "plt.ylabel('g-r')\n",
    "plt.xlim(2.5, 4)\n",
    "\n",
    "plt.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For more on what matching options are available, check out the [separation and matching section of the astropy documentation](http://astropy.readthedocs.org/en/latest/coordinates/matchsep.html).  Or for more on what you can do with `SkyCoord`, see [its API documentation](http://astropy.readthedocs.org/en/latest/api/astropy.coordinates.SkyCoord.html)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercises"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Check that the ``d2d_sdss`` variable matches the on-sky separations you get from comaparing the matched ``coo_sdss`` entries to ``coo_twomass``.  \n",
    "\n",
    "Hint: You'll likely find the ``SkyCoord.separation()`` method useful here."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compute the *physical* separation between two (or more) objects in the catalogs.  You'll need line-of-sight distances, so a reasonable guess might be the distance to HCG 7, which is about 55 Mpc. \n",
    "\n",
    "Hint: you'll want to create new `SkyCoord` objects, but with ``distance`` attributes.  There's also a `SkyCoord` method that should do the rest of the work, but you'll have to poke around to figure out what it is."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Transforming between coordinate systems and planning observations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now lets say something excites you about one of the objects in this catalog, and you want to know if and when you might go about observing it.  `astropy.coordinates` provides tools to enable this, as well."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Introducting frame transformations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To understand the code in this section, it may help to read over the [overview of the astropy coordinates scheme](http://astropy.readthedocs.org/en/latest/coordinates/index.html#overview-of-astropy-coordinates-concepts).  The key bit to understand is that all coordinates in astropy are in particular \"frames\", and we can transform between a specific `SkyCoord` object from one frame to another.  For example, we can transform our previously-defined center of HCG7 from ICRS to Galactic coordinates:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center.galactic"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The above is actually a special \"quick-access\" form which internally does the same as what's in the cell below: uses the `transform_to()` method to convert from one frame to another."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from astropy.coordinates import Galactic\n",
    "hcg7_center.transform_to(Galactic())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that changing frames also changes some of the attributes of the object, but usually in a way that makes sense:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center.galactic.ra  # should fail because galactic coordinates are l/b not RA/Dec"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center.galactic.b"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Using frame transformations to get to AltAz"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To actually do anything with observability we need to convert to a frame local to an on-earth observer.  By far the most common choice is horizontal coordinates, or \"AltAz\" coordinates.  We first need to specify both where and when we want to try to observe."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from astropy.coordinates import EarthLocation\n",
    "from astropy.time import Time\n",
    "\n",
    "observing_location = EarthLocation(lat='31d57.5m', lon='-111d35.8m', height=2096*u.m)  # Kitt Peak, Arizona\n",
    "# If you're using astropy v1.1 or later, you can replace the above with this:\n",
    "#observing_location = EarthLocation.of_site('Kitt Peak')\n",
    "\n",
    "observing_time = Time('2010-12-21 1:00')  # 1am UTC=6pm AZ mountain time"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we use these to create an `AltAz` frame object.  Note that this frame has some other information about the atmosphere, which can be used to correct for atmospheric refraction.  Here we leave that alone, because the default is to ignore this effect (by setting the pressure to 0)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from astropy.coordinates import AltAz\n",
    "\n",
    "aa = AltAz(location=observing_location, obstime=observing_time)\n",
    "aa"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we can just transform our ICRS `SkyCoord` to `AltAz` to get the location in the sky over Kitt Peak at the requested time."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center.transform_to(aa)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Alright, it's up at 6pm, but that's pretty early to be observing.  We could just try various times one at a time to see if the airmass is at a darker time, but we can do better: lets try to create an airmass plot."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# this gives a Time object with an *array* of times\n",
    "delta_hours = np.linspace(0, 6, 100)*u.hour\n",
    "full_night_times = observing_time + delta_hours\n",
    "full_night_aa_frames = AltAz(location=observing_location, obstime=full_night_times)\n",
    "full_night_aa_coos = hcg7_center.transform_to(full_night_aa_frames)\n",
    "\n",
    "plt.plot(delta_hours, full_night_aa_coos.secz)\n",
    "plt.xlabel('Hours from 6pm AZ time')\n",
    "plt.ylabel('Airmass [Sec(z)]')\n",
    "plt.ylim(0.9,3)\n",
    "plt.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Great!  Looks like it's at the lowest airmass in another hour or so (7pm).  But might that might still be twilight... When should we start observing for proper dark skies?  Fortunately, astropy provides a ``get_sun`` function that can be used to check this.  Lets use it to check if we're in 18-degree twilight or not."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from astropy.coordinates import get_sun\n",
    "\n",
    "full_night_sun_coos = get_sun(full_night_times).transform_to(full_night_aa_frames)\n",
    "plt.plot(delta_hours, full_night_sun_coos.alt.deg)\n",
    "plt.axhline(-18, color='k')\n",
    "plt.xlabel('Hours from 6pm AZ time')\n",
    "plt.ylabel('Sun altitude')\n",
    "plt.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Looks like it's just below 18 degrees at 7, so you should be good to go!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercises"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Try to actually compute to some arbitrary precision (rather than eye-balling on a plot) when 18 degree twilight or sunrise/sunset hits on that night."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Try converting the HCG7 coordinates to an equatorial frame at some other equinox a while in the past (like J2000).  Do you see the precession of the equinoxes?\n",
    "\n",
    "Hint: To see a diagram of the supported frames look [here](http://docs.astropy.org/en/stable/coordinates/#module-astropy.coordinates).  One of those will do what you need if you give it the right frame attributes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Wrap-up"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For lots more documentation on the many other features of `astropy.coordinates`, check out [its section of the documentation](http://astropy.readthedocs.org/en/latest/coordinates/index.html).\n",
    "\n",
    "You might also be interested in [the astroplan affiliated package](http://astroplan.readthedocs.org/), which uses the `astropy.coordinates` to do more advanced versions of the tasks in the last section of this tutorial."
   ]
  }
 ],
 "metadata": {
  "astropy-tutorials": {
   "author": "Erik Tollerud <erik.tollerud@gmail.com>",
   "date": "July 2015",
   "description": "Demonstrates use of astropy.coordinates for common tasks. Includes matching catalogs against each other, basic observing planning tasks, and basic usage of coordinates.",
   "link_name": "Using astropy.coordinates to Match Catalogs and Plan Observations",
   "name": "",
   "published": true
  },
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython"
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
